/*
 *  lib-smv: Matrix - Vector Library targeted to Structural Mechanics Problems
 *  Copyright (C) 2006-2008  St.Brcic, Lj.Zugic
 *
 *  This program is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */
 
/** @file TQLI.h
 *  @brief Definition of class TQLI - solution of the eigenvalue problem
 *
 *  EIGENVALUE PROBLEMS:  Class TQLI   
 *  STANDARD EIGENVALUE PROBLEMS A*x = L*x
 *  
 *  Householder transformation to tridiagonal form
 *  QL transformation with implicit shifts
 */ 

#ifndef _TQLI_H
#define _TQLI_H

#include <iostream>
#include <fstream>
#include <iomanip>
#include <cstdlib>
#include <cstdio>
#include <cmath>
#include <cassert>
#include <string>
#include "Arrays.h"
#include "iArrays.h"

using namespace std;

namespace smv_lib
{

//! utilities
#define SIGN(a,b) ((b) >= 0.0 ? fabs(a) : -fabs(a))
#define DSQR(a) (a) == 0.0 ? 0.0 :((a)*(a))
//! computes (a2 + b2)^1/2 without destructive underflow or overflow
double pythag(double a, double b);

/**
* \brief class TQLI
*
*  STANDARD EIGENVALUE ANALYSIS
*   A*x = L*x
*
*  Householder Reduction combined with
*  QL Transformation with implicit shifts
*  (mainly according to "Numerical Recipes in C").
*/
class TQLI
{
  private:
    //! Matrix A for the standard eigenvalue problem A*x = Lambda*x
    Matrix A;
    //! Order of matrix A
    int N;
    //! Modal matrix: columns of V are eigenvectors
    Matrix V;
    //! Vector of eigenvalues: Spectral vector (matrix)
    Vector d;
    //! Indicator: true if matrix A is given - for standard problem
    bool isMatrix_A;
    //! Indicator: true if the problem is solved
    bool isCalculated;
    //! Name of the modal matrix
    string V_Name;
    //! Name of the spectral vector
    string d_Name;
    //! Indicators: true if names are given
    bool isVName, isdName;

  public:
    //! empty (default) constructor
    TQLI(void);
    //! constructor for the standard complete problem
    TQLI(Matrix& a);
    //! destructor
    ~TQLI(void);

    //--------- set and get functions (in case of a default constructor)
    //! set Matrix A
    void setMatrix_A(Matrix& a);
    //! set modal matrix name
    void setV_Name(string vname);
    //! set spectral vector name
    void setd_Name(string dname);
    //--------- calculate eigenvectors and eigenvalues
    void calculate();
    //--------- give the results
    Matrix& giveModalMatrix();
    Vector& giveSpectrum();

  private:
    //! Householder reduction to tridiagonal form
    void tred2(Vector& e);
    //! QL algorithm with implicit shifts (after Householder reduction)
    void tqli(Vector& e, bool isTred2 = true);
    //! sort obtained eigenvecetors / eigenvalues
    void eigsrt_tqli();
};

//=================== external functions
//! Householder's reduction
void tred2  (Matrix& A, Vector& d, Vector& e);
//! QL algorithm with implicit shifts
void tqli   (Vector& d, Vector& e, Matrix& A, bool isTred2 = true);
//! sort the results
void eigsrt_tqli (Vector& d, Matrix& V);

}  // end of namespace smv_lib

#endif //  _TQLI_H

//=============================== end of file
